Response of the Shockley surface state to an external electrical field: 
A density- functional theory study of Cu(lll) 



K. Berland 1 -*, T.L. Einstein 2 ^, an d P. Hyldgaard 1 ^ 

1 Department of Microtechnology and Nanoscience, MC2, 
Chalmers University of Technology, SE-41296 Goteborg, Sweden 
2 Department of Physics, University of Maryland, College Park, Maryland 20742-4111, USA 

(Dated: December 20, 2011) 

The response of the Cu(lll) Shockley surface state to an external electrical field is characterized 
by combining a density-functional theory calculation for a slab geometry with an analysis of the 
Kohn-Sham wavefunctions. Our analysis is facilitated by a decoupling of the Kohn-Sham states via 
a rotation in Hilbert space. We find that the surface state displays isotropic dispersion, quadratic 
until the Fermi wave vector but with a significant quartic contribution beyond. We calculate the 
shift in energetic position and effective mass of the surface state for an electrical field perpendicular 
to the Cu(lll) surface; the response is linear over a broad range of field strengths. We find that 
charge transfer occurs beyond the outermost copper atoms and that accumulation of electrons is 
responsible for a quarter of the screening of the electrical field. This allows us to provide well- 
converged determinations of the field-induced changes in the surface state for a moderate number 
of layers in the slab geometry. 



PACS numbers: 73.20.At,71.15.Mb,73.90.+f 
I. INTRODUCTION 

Surface statesPH are long known to significantly affect 
the properties of surfaces in a host of ways. This is par- 
ticularly so when the surface state crosses the Fermi level, 
rendering it metallic. Such states provide the possibility 
of low-energy adsorption and enhancement of transport. 
On the close-packed (111) surfaces of noble metals, such 
states have their minima at the center of the surface Bril- 
louin zone (r) in the bulk L gap, have minimal angular 
anisotropy, and are well approximated by free-electron 
dispersion!^ Since the energy difference between the 
Fermi level and the bottom of the band is small, so is the 
Fermi wave vector, leading to a Fermi wavelength that is 
nearly an order of magnitude larger than that of typical 
bulk Fermi wavelengths. Furthermore, confinement to 
the surface leads only to slow decay in directions perpen- 
dicular to the surface plane, Fig. [T] Arguably the most 
dramatic outcome is the observation of quantum corrals 
and mirages on surfacesP^ These states can also pro- 
duce ordered superstructures of adsorbed species P^ELl] 

A question of both fundamental and practical con- 
cern is the sensitivity of these metallic surface states to 
perturbations. Moderately strong perturbations such as 
chemisorption can easily destroy the surface statep2322] or 
create overlayer resonances pDESl There is also a theoret- 
ical prediction that alkali-overlayer formation can lead 
to a localization (for dynamics parallel to the surface) 
of high-energy electrons in a resonance stated On the 
other hand, weak perturbations can allow the state to 
survive but with an altered dispersion relation (typically 
a shifted minimum and a changed curvature) 121ES1 This 
offers the exciting possibility to manipulate the Fermi 
wavelength and effective mass of the state, if one can un- 
derstand in detail how perturbations such as adsorption 
affect dispersion. 




FIG. 1: (Color online) Contours of the F-point surface- 
state density. Yellow (red, blue) surfaces correspond to high 
(medium, low) density. The density is largest outside the sur- 
face, and lowest in the planes formed by the atoms. Graphics 
generated with VESTAP 



However, adsorption typically leads to several differ- 
ent changes to the surface.^ First, there can be charge 
transfer, producing an electric field due to the resulting 
surface dipole. (There could also be effects from the in- 
trinsic dipole of an organic adsorbate.) Second, there 
can be correlated electron hopping such as characteris- 
tic of covalent bonds. Third, the adsorbate will perturb 
the tails into the vacuum of the metal-surface electron 
density.^ Furthermore, as dipole-producing adsorbates 
approach each other, they create a depolarizing field that 
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decreases the dipole, consistent with experiments for Na 
atoms on Cu(001)pS Hence, the prospect of predicting 
how a particular adsorbate modifies the surface state 
poses a considerable challenge. 

In this paper, we instead focus on just the first as- 
pect of the adsorption bond, that due to charge transfer. 
Specifically, we look at the effect that a uniform electric 
field normal to the surface has on the surface state. This 
problem also has advantages over a direct study of a pure 
ionic bond: the electric field is uniform, and the screen- 
ing is well confined to the direction perpendicular to the 
substrate. Furthermore, this simple scenario bears on 
the more complicated geometry involved with th e effe ct 
of strong fields between STM tips and substrates! 29 * 30 ! 

For over four decadesp^with increasing sophistication, 
theoreticians have used density- functional theory (DFT) 
in various implementation^S"^ ^ examine the effect of 
static perpendicular electric fields on the electronic prop- 
erties of simple metal surfaces (until relatively recently, 
typically jellium-like). Trying to understand experiments 
related to self-diffusion on Pt(OOl), Feibelman discussed 
how such electric fields can be used to tune diffusion bar- 
riers and possibly alter the dominant mechanism of mass 
transport.^ Negulyaev et al. showed that such fields 
could serve as a switching tool for magnetic states in 
atomic-scale nanostructuresP^ However, we know of no 
DFT study of the effect perpendicular electric fields on 
metallic surface states. 

We further note that understanding the surface re- 
sponse to a perpendicular electric field and the formation 
of an image plane are fundamental building blocks in the 
study of the van der Waals interactions. ^"^! The image 
plane is available from DFT calculation d 35 * 47 " 49 ! and the 
handling of screening is a central element in the vdW- 
DF method JSSHSU Knowledge of the image plan e posi tion 
permits an approximation of multipole effects anc j 
hence a quantitative study of interactions of, for exam- 
ple, atoms and molecules at noble-metal surfacesP^^ 
The image plane is also important for ensuring trans- 
ferability over a range of different binding distances, for 
example, in molecular crystals psSHSQl 

Much of our exploration was motivated by a desire 
to manipulate the surface states that mediate the in- 
teractions between anthraquinone (AQ) molecules on 
Cu(lll)P^! the giant regular honeycomb formed sponta- 
neously by them are likely related to such interactions.^! 
Alternatively, the pores in the honeycomb can be viewed 
as an array of two-dimensional (2D) quantum dots. The 
stabization of the pattern may be due to the population — 
from the metallic surface state — of the 2D orbitals of the 
dots, forming what amounts to closed-shell, 2D-noble- 
gas-like quasiatoms!^ 3 -! By manipulating the surface state, 
we hope to tune its Fermi wave vector (and, concomi- 
tantly, its effective mass) and thereby enhance or desta- 
bilize different superlattice structures. A major goal of 
this study of surface-state response is to gain the ability 
to engineer novel structures on surfaces. Furthermore, 
the standing waves within the honeycomb cells, arising 



from these surface states, are believed to determine the 
potentials that small molecules like CO encounter when 
adsorbing within the cellsP^l These cells form a set of 
identical nanostructures with thermodynamic-like behav- 
ior that differs significantly from that found when these 
molecules adsorb on large defect-free flat surfaces. 

Most traditional DFT implementations rely on super- 
cells to model surfaces, with slabs separated by vacuum. 
However, since pairs of Shockley surface states on oppo- 
site sides of a slab hybridize, the energy and dispersion 
are affected, with ensuing loss of accuracy. Use of very 
thick slabs can marginalize the hybridization effects, but 
such a brute-force approach is computationally expen- 
sive and more susceptible to numerical noise. In this 
paper, we present a method which extracts proper, un- 
hybridized, surface states from standard supercell-DFT 
studies for geometries with moderate slab thicknesses. 
Our method is based on a simple rotation in the Hilbert 
space spanned by the two Kohn-Sham (KS) metallic sur- 
face states that are found in underlying semilocal DFT 
calculations. We use our method to characterize the re- 
sponse of the Cu(lll) surface state, but it should also 
be useful for the study of other, more complex, material 
systems. 

The KS-rotation method presented here is complemen- 
tary to the use of more advanced DFT implementations, 
such as the embedding Green function methodf^HH! 
which effectively model a semi-infinite surface. The em- 
bedding method can also handle an external fieldP^l As 
our study is based on a traditional DFT implementa- 
tion, we do not correctly describe t he im age-potential 
behavior, for which GW calculation; ] 72 ! 73 ! are normally 
required. The embedded method allows an explicit 
inclusiorP ^ 68 ! 71 ! of an image- like behavior and can there- 
fore determine image-potential states in DFT. The inclu- 
sion of this image-like behavior would likely also improve 
the accuracy of the evanescent part of the surface states 
found in the slab analysis. In spite of these benefits of 
an embedding method, we believe it is also important 
to continue to seek simple mechanisms to enhance the 
accuracy in widely used (for examples, Refs. I23I74H55]) 
slab-geometry DFT studies of surface states. 

The plan of this paper is as follows. In Sect. II we dis- 
cuss our computational methods, paying particular at- 
tention to a simple yet unambiguous and robust way to 
accurately decouple the surface states on the two sides 
of the slab and to obtaining a faster convergence as a 
function of slab thickness. In Sect. Ill we character- 
ize the surface state with no electric field, while in Sect. 
IV we describe the changes in the wavef unctions, poten- 
tial profile, and dispersion in the presence of an applied 
perpendicular field. Sect. V discusses the screening of 
this electric field, deriving the relative contribution of 
the electrons in the surface state. It also makes compar- 
isons with experimental data and illustrates the decou- 
pling method for benzene on Cu(lll). Finally, Sect. VI 
offers conclusions about the impact of these changes in 
the surface state. 
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II. COMPUTATIONAL METHODS 

The electronic structure is obtained with DFT 
within the generalized-gradient approximation (GGA) 
for exchange-correlation using the PBEP version. For 
these calculations, we use the ultrasoft pseudopotential 
plane- wave code DACAP0,^ ! with an energy cutoff of 400 
e ylH] anc j a k-sampling of 16x16x1. The KS states re- 
quired to obtain the surface-state dispersion is calculated 
non-self consistently in a post-processing Harris func- 
tional calculation 9 ^ in which the wave vectors were sam- 
pled with a grid spacing of 0.01 K, where K is the size 
of the shortest in-plane reciprocal-lattice vector. 

The Cu(lll) surface is modeled as a finite slab in a 
supercell since the plane-wave scheme restricts us to us- 
ing periodic boundary conditions. The top and bottom 
copper layers in two adjacent supercells are separated by 
12A, thus insuring negligible cross-coupling. The lattice 
constant of the slab is set to that of copper, a = 3.65A, 
as obtained in a separate bulk PBE calculation. The 
electrical fields used in this study only slightly perturb 
the electronic structure of the surface, and we find no 
significant relaxations of the atoms of the surface slab; 
the atoms are therefore frozen in their truncated bulk 
positions! 9 ^ A dipole layer in the vacuum region in duces 
an external electrical field normal to the copper slab! 94 * 95 -! 
Both positive and negative electrical fields are studied in 
a single calculation. 

A. Issues with coupled surface states 

The finite slab geometry makes the surface state couple 
to both the two sides of the slab, and therefore challenges 
the analysis of surface state response to an external per- 
turbations. 

In the upper panel of Fig. [2j the full curves show 
the sum-projected density of the KS wavefunctions with 
surface-state character, ipf^ s (k) at the T point, 

pf s (z) = J Ax J &y\^f s {x,y,z)\\ (1) 

for a six-layer thick copper slab in zero external field. 
These states couple equally to both sides of the slab, 
rather than being localized on one side. They therefore 
lack the characteristic exponential decay into the bulk. 
In a hybridization, or tight-binding picture, these KS 
states can be viewed as linear combinations of surface- 
localized (SL) states that hybridize in a finite slab geom- 
etry. For zero electrical field, the KS states form sym- 
metric and anti-symmetric combinations of the underly- 
ing SL states (which provide an accurate description of 
the actual surface-state behavior). 

The dashed lines in the top panel of Fig. [2] shows the 
sum-projected densities, pf L , of these SL states. The SL 
states exhibit "good" surface-state properties, like expo- 
nential decay into the bulk and localization at the sur- 
face. For sufficiently large slabs, these SL states are good 



representations of the proper surface states, i.e., surface 
states as they would have been calculated in an accurate 
DFT of a semi-infinite bulk system. 




£[V/A] 

FIG. 2: (Color online) Issues with surface-state hybridization 
that arise because of the finite slab geometry. The upper panel 
shows the sum-projected surface-state densities obtained for 
a six layer wide copper slab. The full lines shows the Kohn- 
Sham (KS) states with surface-state character. The dashed 
lines show the surface-localized (SL) states for which we here 
present a simple but robust determination (since the SL states 
provide a more accurate description of genuine surface-state 
behavior). The lower panel shows the electric field varia- 
tion of the surface state for a 15-layer slab calculation using 
the plane-wave code DACAPO. The full line shows the least- 
squares fit to decoupled surface-state energies as obtained for 
the SL states, Sec. IV.B. 

In the lower panel of Fig.[2]the filled (red) circles show, 
for different external fields, the KS eigenvalues for the 
states with surface-state character, ef- s . These results 
were obtained for a 15 layer thick slab. The jump at 
E = in the KS calculation of the minimum surface- 
state energy, e (0), indicates an avoided crossing, which 
further supports the picture of a coupled two-level sys- 
tem. Since the field influences both the properties of the 
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underlying states and the linear combination making up 
the KS states, it is a challenge to deduce the inherent 
response of the proper surface states. 

B. Decoupling of surface states 

This subsection presents our numerically robust 
method to construct the (pair of) SL states from the 
KS surface states of the slab. The full line in the bot- 
tom panel of Fig. 2 shows the effective Fermi level shifts 
for the SL result. The continuity of this result contrasts 
the singular response which appears in an analysis based 
directly on the KS states (shown as [red] dots). 

Our assumption is that the KS surface states arise ex- 
clusively from hybridization of the actual SL states; we 
assume that these do not couple to bulk states or surface 
resonances. In this case, we can for a given k set up a 
bonding/antibonding Hamiltonian for this SL two-level 
system: 

(e(k) + w(k) n(k) \ rl 

H{K)-y e (k)-W{k)J [> 

Here e(k) ± W(k) gives the energy of the uncoupled 
states, where W is the detuning between the levels and 
f2 is the coupling. In terms of them, the eigenvalues of 
the KS wavefunctions are given by 

<^ = e±y/W*+W, (3) 

which follows from the invariance of TrH and DetH. As 
the slab size increases, fl — > and the eigenvalues of the 
decoupled surface states reduce to e ± W. 

The key observation is that the SL state energy e con- 
verge much faster as the slab size increases than the cou- 
pling parameter f2 vanishes. The same is generally true 
of the detuning parameter W; even at E — (when the 
correct detuning must vanish identically), we find (below) 
that the numerical value of W converges fairly well before 
the accumulation of numerical noise eventually destroys 
this behavior. Therefore, by constructing T-L of Eq. ([2]), 
we can characterize the surface-state properties to high 
accuracy, with a smaller unit cell than what is required 
for the KS to decouple due to asymmetry induced by 
small-scale numerical effects. 

In general a SU(2) rotation is needed to construct the 
SL basis states and T-L as we seek a linear combination of 
the KS wave- functions ipi(k) that localizes the state on 
one side of the slab through constructive and destructive 
interference. However, when we can identify the symme- 
try point of the unit cell in the xy-plane, we can ensure 
that both KS states have the same variation in the com- 
plex plane. This is done, for example, in the present 
Cu(lll) slab study, by setting the KS states real at this 
symmetry point via multiplication by a simple phase fac- 
tor: ipi — > exp(ir])ipi. The transformed KS states must 
then have the same complex-phase variation across the 
unit cell because the KS states are Bloch functions: they 



can be written as ipi(r) — i*i(r)e lkr , where Ui(r) has 
the same periodicity as the supercell and can be chosen 
realP 

With this transformation, the full SU(2) transforma- 
tion reduces to an 0(2) rotation R in the Hilbert space: 

( <t>i \ _ ( cos0 -sin 6» \ ( ip x \ , . 

\<h ) \ sine cose J \ yj 2 J ■ w 

To determine the value of 6*, we choose a condition of 
maximal localization (MaxLoc), that is, we minimize the 
sum of the variances of the generated wavefunctions 

F(e) = Y,m\^)\-m\z*m 2 - (5) 

i 

Such a condition has also been used to construct the 
generalized- Wannier functions.^ 

The matrix transformation that corresponds to this 
rotation is H -> RHR T = H S l- Thus, 

f ef s \ /e KS + iAe KS cos20 iAe KS sm20 \ 
V 4 S / V §Ae KS sin20 e KS - \ Ae KS cos 20 ) ' 

(6) 

Here e KS = (ef s + dp)/2 and Ae KS = (ef s - ef s )/2 By 
comparing to Eq. Oj), we identify the SL-hybridization 
parameters, e = e , W = (Ae K s/2) cos(20), and ft — 
(Aexs/2) sin(20). We note that fl is real because we can 
here work with an 0(2) rotation. 

For perfectly symmetric surfaces we have e = 7r/4, 
and the surface-state energy is then the average of the 
eigenvalues of the two states. Thus, the full machinery 
discussed here is not required. For certain perturbations, 
like adsorbate systems, a less elegant solution is to use 
symmetric adsorbates on both sides of the slab. Such 
a brute-force approach has the drawbacks of increased 
computational costs and fewer layers with a bulk-like be- 
havior; furthermore, wave functions are not decoupled. 
Strong asymmetric perturbations, such as halogen over- 
layers, lead to W ^> Q and a natural decoupling of the 
states. However, a natural decoupling does not happen 
for weak perturbations, like adsorbates bound by van der 
Waals interactions,^ or for dilute chemisorbed overlay- 
ers, where £1 ~ W unless a huge number of layers are 
used, typically beyond the computational feasibility. 

C. Convergence of surface-state properties 

Our method to construct the SL states significantly 
reduces the number of layers needed to accurately char- 
acterize the surface state; however, the slab still must 
be thick enough to describe most of the relatively slow 
decay into the bulk. We also confirm that the MaxLoc 
condition (see Eq. [5| properly decouples the two states 
by comparing with the results for a 24 layer calculation, 
where decoupling arises numerically as Q — > 0. 

Table |T] shows the calculated SL state parameters for 
for different slab thicknesses at zero electrical field. The 
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TABLE I: Convergence of the surface state with slab thick- 
ness. All three columns of energies are in meV. The bold 
numbers, for 15 layers, are used for the rest of this paper. 



# layers 


e F - e(0) 




W (E=0) 


6 


520 


272 


0.24 


9 


464 


113 


0.38 


12 


449 


50 


0.17 


15 


443 


23 


0.81 


18 


443 


10 


1.34 


24 


442 


2 


11.0 



Fermi energy ep relative to e(0) converges to the sub- 
meV level for a 15-layer slab if the decoupling method 
is used. For a six-layer slab the surface-state energy dif- 
fers from the converged value by 80 meV, or about 20%; 
we consider it a minimum slab thickness for an approx- 
imate account of the surface state, which can be useful 
for studying surface-state shifts for adsorbates-systems 
requiring a large supercell in the in-plane direction. 

We find that coupling SI decays significantly slower 
than the value of ep — e(0) converges, a fact which as 
mentioned above motivates our approach. We note that 
the value of 272 meV for six layers will almost deplete 
one of the two surface-related KS bands (while driving 
the other at least partially into the energy range of bulk 
states). Nevertheless, the corresponding result for the 
effective depth of the surface-state, Fermi see, ep — e(0) 
differs from the converged value by merely 57 meV. For 
15 layers fi reduces to 23 meV, while the cf — e(0) is 
converged to within less than an meV. 



D. Effects of numerical noise on optimal slab 
geometry 

The nonzero value of the detuning W at zero electri- 
cal field stems purely from numerical noise and grid ef- 
fects, since the slab geometry is symmetric. The upper 
panel of Fig. [4] shows the KS states for different electri- 
cal fields. For zero electrical field, the slight asymmetry 
of the curves hints of the nonzero W. For larger fields, 
the KS states clearly favor one side of the slab, but even 
for E = 0.36V/A, the wavefunctions are localized on 
both sides of the slab. The inserts illustrate (purple dot) 
the coordinate in the 2D space (W,fl), with magnitude 
(e^ — e)p)/2, constructed using the MaxLoc condition. 
The detuning grows as the electrical field increases, while 
the coupling remains roughly constant. 

For six layers the two KS surface-state densities are 
fully symmetric as shown in the upper panel of Fig.[2j For 
even wider slabs, the coupling reduces, and the surface- 
state densities becomes increasingly asymmetric; for 24 
layers, we find that the surface states localize almost ex- 
clusively on one side of the slab, consistent with a detun- 
ing five times larger than the coupling SI. 



That e s « e for 24 layers reflects this finding and 
the general expectation that increasing the number of 
layers must eventually decouple the KS and ensure an 
automatic approach to the more meaningful SL states. 
Nevertheless, we find that the increasing value of W for 
slabs with than 15 layers is the sole reason for this de- 
coupling. As this increase in W (for E = 0) is solely an 
expression of an accumulation of numerical noise (with 
the number of atoms and electrons), we conclude that 
such brute-force decoupling is not desirable. 

We propose instead to use the maximum-localization 
approach to identify the genuine surface-state behavior 
and avoid such numerical noise. In this study we use 15 
layers, as it is an optimum between a converged value of 
e and a minimal numerical noise; our approach converges 
the surface-state energies to meV accuracy. 

III. SURFACE STATE DISPERSION IN ZERO 
ELECTRICAL FIELD 

Our calculations of the surface-state dispersion in zero 
electrical field are pre sented in this section and com pared 
to experimenta l 99 * 10 -^ and earlier theoretical studiesP 5 * 101 ! 
We find that the surface state is practically isotropic even 
for k = |k| > 2.5kp, and that non-parabolicity becomes 
significant once k exceeds kp. 

Figure [3] shows the surface-state dispersion as a func- 
tion of the absolute value of k. The filled (green) circles 
indicate the calculated values of e(k) — ep, which have 
been sampled evenly on a k x , fc y -grid, while the (red) 
crosses indicate the energies obtained for k y = 0. That 
they align to form a curve shows that the surface state is 
isotropic, even as far as ~ 2.5|kp|. This isotropy is fur- 
ther evidenced by the circular constant-energy contours 
displayed in the insert. The full curve gives the parabolic 
dispersion, while the dashed includes non-parabolicity 
via a quartic term: 

<k) - e(0) = h ^-ak i + 0(k 6 ) , (7) 
Ira 

with parameters obtained as described in the following 
paragraph. These curves show that up to about the kp, 
the dispersion is well described by a parabolic form, but 
for larger wave vectors, non-parabolicity is significant. 

The deviation from a parabolic-dispersion behavior has 
been observed experimentally and can b e un derstood in 
terms of an s-band tight-binding modelP^l We expand 
e(k) as given by the standard Hamiltonian and find: 

e(k) - e(0) cx 1 - i cos(k x a) - ? cos {^^j cos 
1/, n9 1 /, ,i 10 + cos(66>) „ Nfi 

(8) 

where 9 = arctan(/cj,/fc a .) and a — 2.58A is the nearest- 
neighbor distance. Thus, the quartic correction is neg- 
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ative (albeit much smaller than found in the DFT cal- 
culation), and anisotropy does not appear until the fc 6 
termP^H 

Table [XT] gives deduced surface-state properties and 
compares them to other studies. In obtaining these val- 
ues, we only needed data points obtained for k y = 0, since 
we have amply demonstrated the isotropic surface-state 
dispersion; this procedure also avoids excessive weight- 
ing of large- |k| values, as the number of 2D grid points 
grows approximately linearly with the wave-vector mag- 
nitude. Our values for the effective Fermi level is similar 
to earlier experimental and calculated values. Our value 
is somewhat smaller than that given in Ref ITU which is 
based on a semi-infinite approach including an image po- 
tential. Our value for the effective mass is smaller than 
those of the earlier studies, which is partly a result of the 
fitting procedure used to obtain its value. 

A parabolic fit is often used to extract the effective 
mass from the dispersion curve; however, the effective 
mass is defined by the second-order Taylor expansion of 
the dispersion curve, not by an optimal parabolic fit to 
a curve that may deviate from parabolicity. To properly 
extract the terms in the Taylor expansion with a polyno- 
mial fit, we rely only on values fairly close to the Y point 
and, in addition, include higher-order terms in the fit to 
minimize their influence on the extraction of the lower- 
order ones. It is possible to avoid the influence of quar- 
tic terms on the effective mass with a purely parabolic 
fit, but this requires that we restrict the domain to data 
points very close to the Y point. This has the disadvan- 
tage that the influence of noise is larger. In detail, our 
approach is as follows: first, the effective mass m is ob- 
tained with least-squares using a fourth-order polynomial 
fit for the six smallest k points (corresponding to 6% of 
the reciprocal vector); next, we keep the mass fixed and 
determine a using a sixth-order polynomial using all 20 
data points. 

The correspondence with the data points of Fig. [3] for 
small and medium |k| corroborates our procedure, as just 
described. If we instead include ten data points and fit to 
a purely parabolic dispersion, we find a mass of 0.38m e , 
closer to that of earlier studies (listed in Table |Ti]) . Thus, 
that our deduced mass is smaller than previously ob- 
tained values relates to our use only of data points close 
to the r point and our inclusion of higher-order polyno- 
mial terms. The quartic prefactor a compensates some- 
what for the smaller mass, and the Fermi wavelength 
Xf and wave vector kp are in good agreement with ear- 
lier results. That sensitivity to fitting domain (parabolic 
fitting) was discussed in Ref [71] They note that when 
comparing with experimental data, it is important that 
same procedure is used in both cases. We argue that, 
ideally, all parameters should be defined in terms of the 
Taylor expansion. When they use a minimal sampling 
of k points around the Y point, they obtain an effective 
mass of 0.303m e , which is closer to, and even smaller 
than the mass we obtain. The importance of making 
sure that higher-order terms do influence the extraction 



of the effective mass is also reflected in the significant 
non-parabolic dispersion that Becker et alP^ found for 
Ag(lll). They extracted the effective mass by consider- 
ing only data point collected close to the f point. 



TABLE II: Properties of the surface state. The value ep — e(0) 
is the difference between the Fermi surface and the min- 
imum of the surface state e(0). Also listed are compar- 
isons of our values of the effective mass, surface-state Fermi 
wavelength and wavevectors, m, Xp,kp. Our DFT value 
for the quartic component of the in-plane surface state dis- 
persion, a, is compared with the (m-dependent) estimates 
(lia) 2 /32m w 1.58(m e /m)eVA 4 given by the expansion ^ of 
the tight-binding behavior observed in Ref. 11021 
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FIG. 3: (Color online) Dispersion of the surface state plot- 
ted as a function of |k| for k sampled evenly on a 2D grid. 
The full (dashed) curve gives the best-fit second- (fourth-) or- 
der polynomial. The dotted line indicates the Fermi surface. 
The filled [green] circles, which form a thick line, indicate the 
calculated values, while the [red] crosses indicate those for 
k y = 0. The insert shows the energy contours as function of 
k, with ticks having the same spacing as in the main figure. 
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IV. SURFACE STATE IN AN EXTERNAL 
ELECTRICAL FIELD 

The dispersion of the Shockley surface state was char- 
acterized in the previous section at zero electrical field. 
In this section, we present the results for a finite external 
electrical field. 



accumulating at the dipole layer, ruining the physical 
picture of semi-stable surface states. 



A. Wavefunctions and potential profile 

Figure [T] gives intersections or contours of the variation 
of the surface-state density in three dimensions. Fig. [T] 
confirms that the surface state is mostly located between 
the copper layers or outside the outermost copper layer. 
The figure also shows that, in the in-plane direction, the 
density is largest close to the copper atoms. 

The upper panel of Fig. [3] shows the sum-projected 
density of the KS wavefunctions at the T point, 
(z), for different external electrical fields E. For 
zero field (top line), the two wavefunctions form approx- 
imately symmetric and anti-symmetric functions as de- 
picted by the full (green) and dashed (black) curve, re- 
spectively. That the two curves have almost equal mag- 
nitude on the left side, but not on the right, reflects the 
presence of numerical noise and grid sensitivities in the 
plane-wave DFT calculations. 

The lower panel of Fig.Hlshows the sum-projected den- 
sity of the SL states at the T point, p i (z), and the 
average potential energy as a function of z. The SL states 
differ marginally for different electrical fields; hence, we 
only display those for zero electrical field. These states 
exhibit oscillatory exponential decay into the bulk and 
are located solely on one side (or the other) of the slab, 
further validating use of the MaxLoc condition. The filled 
(yellow) circles indicate the location of the copper atoms. 
The first node of the wavefunctions coincides with the po- 
sition of the top copper layer. The (laterally averaged) 
potential energy profiles for the three field strengths are 
given by surfaces of semi-transparent grey shading; thus, 
the light-, medium-, and dark grey areas indicate ener- 
gies which are smaller than one, two, or all three of the 
potentials. A sharp edge between dark grey and white 
indicates that the potentials overlap, for example within 
the slab. At the left side of the slab, the finite electri- 
cal fields go outward from the surface, corresponding to 
a total reduction of electrons. The number of surface- 
state electrons are also reduced on the left side as the 
electrical field push the electrons toward the surface and 
thereby increases confinement. This raises the minimum 
energy e(0) of the surface state, lowering ep — e(0). At 
the right side, electrons accumulate and so do the surface 
state electrons; they experience weakening confinement, 
hence a lowering of the energy which leads to an increas- 
ing value of e_F — e(0). 

For large electrical fields, the surface state becomes in- 
creasingly unstable, and at some point, depending on the 
size of the vacuum region in the unit cell, electrons start 




FIG. 4: (Color online) The average potential profile and KS 
and SL surface-state densities at the f point for a Cu(lll) 
surface in an external electrical field E pointing to the left, 
i.e., inwards to the surface at the top (at z ~ 0) of the slab 
and outwards from the surface at the bottom (at z ~ — 30A). 
Thus, the force on an electron is away from the surface at 
the top of the slab and the minimum energy e(0) of this top- 
surface SL state is lowered, increasing €f — e(0). The upper 
panel shows the KS surface-state densities (full and dashed) 
for three different strengths of E. At zero electrical field they 
are almost symmetric, while for larger field they almost de- 
couple. The inserts illustrate the corresponding hybridization 
parameters. The lower panel shows the atom positions, filled 
(yellow) circles, and the decoupled SL surface states, evanes- 
cent oscillatory curves. The bottom panel also shows, by sur- 
faces of semi-transparent grey shading, the potential profile 
for the three different strengths of electrical fields. Field- 
induced changes in the potential profile are observable only 
outside the outermost Cu atoms. 
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Figure [5] shows the calculated values of the surface- 
state energy relative to the effective Fermi level (upper 
panel), the effective mass (mid panel), and in the lower 
panel, the relative shift in the wavelength and the wave 
vector based on the shift in Fermi level and effective 
mass in the parabolic free-electron gas approximation. 
It shows that for the largest plotted electrical field the 
wavelength shifts by about 4%. For the upper and mid 
panel, the line gives the least-square fit according to lin- 
ear relations: 



- <0)] E 



[(F - e(0)] E=0 + 
^/m e + A fn E 



A f E 



m 



(9) 
(10) 



where rh denotes m/m e . 

The parameters determined are listed in Table |III| 
along with the shift in the characteristic Fermi wave vec- 
tor fcp (Ak F ) and in the wavelength Xp (A\ F ). 



TABLE III: Response of the surface state to the external elec- 
trical field, obtained from fitting to the calculated dispersion. 
The values A e , A^, A\ F , and Ak F are our results for the field 
derivatives of the surface-state Fermi sea level, £f — e(0), effec- 
tive mass, Fermi wavelength, and Fermi wavevector, respec- 
tively. The value of A a is our DFT-based determination of 
the fraction of overall external-field screening which originates 
from the surface-state behavior. 



A e [e A] 
A A \k/V] 

A Aj ,[AV 

Aa [(VA) _ 



-0.146 (±0.007) 

0.029 
-0.025 

3.68 

-4.82 • 10~ 3 



FIG. 5: (Color online) Calculated energy at the F-point (up- 
per panel) and mass (middle panel) shift for different exter- 
nal electrical fields. Here we take a positive value of the field 
strength E to imply that the field E (the force) points away 
from (towards) the surface. In terms of the underlying slab 
calculations, Fig.|4j the positive- E (negative- E) results char- 
acterize the behavior of the SL state at the left-most (right- 
most) slab surface. The filled (red) circles [(white) diamond] 
gives the calculated values for a 15 [12] layer slab. The black 
line gives the least-square fit to the calculated data. In the 
lower panel the relative shift in wavelength (wave vector) is 
given by the full (dashed) curve. 



B. Modified dispersion 

The external electrical field not only modifies the effec- 
tive Fermi level €f of the surface state, but also influences 
the surface-state dispersion, expressed here in terms of a 
shift in the effective mass. The mass is obtained as de- 
scribed in the previous section. 



The range of field strengths E plotted in Fig. [5] are 
larger than values which are directly achievable in most 
experiments (but smaller than what is regularly achieved 
in adsorption studies, as discussed in Sec. V). In partic- 
ular, in a study of the Stark shift of the surface state 
due to measurement by scanning tips rather than pho- 
toemission, Berndt's grourj^l found a linear decrease in 
[ £ f - e(0)] B vs. E. This is true up to E = 0.055 V/A, 
at which there is a energy shift of magnitude ^13 mcV; 
thereafter, the rate of change increased markedly, an ef- 
fect which the authors attributed to the breakdown of 
the tunneling regime at small tip-surface separation. We 
note that while the experimentally observed breakdown 
occurs around 0.055 V/A for Cu(lll), it occurs already 
at 0.008 V/A for the more weakly bound Shockley state 
on Ag(lll). 

From their graph of energy shift vs. field, we extract 
the mean slope to be A t ss —0.23 eA for Cu(lll), com- 
parable to our computed value in Table III We speculate 
that the larger magnitude in the experiment arises from 
the additional strain in lateral directions due to the non- 
uniform field between the STM tip and the surface. 
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V. DISCUSSION 

A. Role of the surface state in screening an 
external field 

In a metal, the electrons close to the Fermi level com- 
pletely screen an external static electrical field perpen- 
dicular to the surface. The virtually identical potential 
profile inside the slab, for — 30A < z < 0, in the lower 
panel of Fig. [4] for different electrical fields shows that 
our DFT calculations accounts for this screening. 

Figure [6] displays the charge density response to the 
external fields. The charge is induced outside the outer- 
most copper atom, while there are small oscillations into 
the first two bulk layers. The response is almost iden- 
tical, with opposite prefactor, at the other side of the 
slab. The short screening length indicates that the bulk 
electrons play a major role in screening the electric field. 




FIG. 6: (Color online) Charge transfer induced by external 
electrical field averaged over the in-plane directions. The 
upper light (cyan) [dark (black), dashed (orange)] is for 
E= -0.36 [-0.12, 0.24] V/A. 

To determine the fraction of the screening performed 
by the surface state compared to the bulk states, we ob- 
tain the total charge induced on the surface from Gauss's 
law and the charge accumulated in the surface state with 
its parabolic free-electron dispersion. Gauss's law re- 
lates the induced charge on the surface to the electri- 
cal field by o — e^E. To obtain the charge accumula- 
tion in the surface state, we recall that 2d-electron den- 
sity of states (including spin degeneracy) is mtp /ir (as- 
suming in this paragraph, for notational simplicity, that 
e(0) = 0). Then the charge density due to the surface 
state is a s — —emep/ir at T = 0, which is the appro- 
priate temperature for comparison with standard DFT 
calculations. Since the mass m and the Fermi-level €f 
change as an external electrical field is applied, so does 
the charge of the surface state. To first order in the elec- 
trical field, this shift is given by A<r s = —A(meF)e/^ = 



-A^eE/n + 0(E 2 ), where A a = (mA t + Arae F ). Com- 
bining the expressions, we find the fraction of the elec- 
trical field screened by the surface state: 

— = ^^+(D(E) = 0.27. (11) 

The surface state therefore plays a significant but not 
dominant role in the screening of the electrical field. We 
note that the calculated potential profile and induced 
charge density curves exhibit minute oscillatory variation 
for different fields. This observation, again, accords well 
with the idea that bulk electrons perform most of the 
screening. 

That screening arises predominantly from bulk elec- 
trons has several noteworthy consequences. It sup- 
ports the practice of calculating interactio n energ ies due 
to surface-state mediated interaction jlUtl21 15 l 17 l b ase d 
solely on the changes in the energies of single-particle 
states, ignoring the electrostatic term in the full Harris 
formalism!^ It also suggests that interactions between 
closely spaced adsorbates (within a few lattice spacings 
of each other) are dominated by the bulk states rather 
than the surface states. Based on a STM study Petersen 
et alP^S concluded that screening of step edges at the 
surface was dominated by the surface states on Cu(lll) 
an Au(lll), but they also noted that screening of defects 
slightly below the surface is dominated by bulk electrons. 
The important role of the bulk electrons is also reflected 
in the li near di spersion of acoustic-surface plasmons on 

cu(m)P!iina 

B. Comparison with experiments 

In Table |IV| we collect results for adsorption-induced 
shifts in tp — e(0), as well as changes in the effective 
mass and the workfunction when available, for Cu(lll) 
and Ag(lll) to provide information about the kinds of 
values measured in mostly-recent experiments. We have 
not included gold-surface-adsorbate systems because of 
the strong spin- orbit co upling gives rise to a Bychkov- 
Rashba splittingJ^HHI] j n addition, Au(lll), in contrast 
to Cu, Ag, and other fee (111) surface s, reconstructs, tak- 
ing on a herringbone pattern!- 120 * 121 ! We note that is not 
always possible to fully calibrate the values in the exper- 
iments since precise coverages are rarely given. Data is 
typically for a monolayer (1 ML), which refers to the sat- 
uration coverage rather than one adsorbate per substrate 
atom. Even if this information were available, there are 
many other factors, discussed at the outset, which can 
contribute to the shift. This is true especially for non- 
alkali adsorbates. 

Alkali adsorbates invariably lower the surface band, 
eventually dragging it into the bulk continuum, where it 
becomes a resonance; often the details are not reported 
(for example for Cs/Cu(lll)p2l); associated calculations 
are problematic due to the large unit cells needed for 
fractional coverage and ill-defined order. There have 
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TABLE IV: Experimental results for alkali adsorbates and 
electric fields on the Shockley surface state on Cu and 
Ag (111). The ratio of effective to bare electron mass is in- 
dicated by m. Unless indicated otherwise, coverages are a 
saturated overlayer ("1 ML"). The entries for electric fields, 
denoted by "Eisubstrate" , are those discussed in Subsection 
|IVB| the "adsorbed" values are for the the maximum fields 
before the tunneling regime breaks down (not the tabulated^ 
values for R = 500MQ). 



System 


[eF — & 


'0)1 fmeVl 

v. /J L J 




in 




clean 


adsorbed 


clean 


adsorbed 


< §ML Li/Cu °" h 




755 






0.11-0.15ML Na/Cu c 


390 


800 






0.4ML K/Cu d 


410 


755 


0.41 


0.36 


Cs/Cu a 










Ar/Cu e 


434(2) 


376(3) 


0.43(1) 


0.46(3) 


Kr/Cu e 


434(2) 


358(2) 


0.43(1) 


0.44(2) 


Xe/Cu c 


434(2) 


291(2) 


0.43(1) 


0.44(2) 


Xe/Cu / 


440(10) 


310(10+) 


0.40(2) 


0.42(3) 


\ML Na/Cu 9 


~340 


690 






C 6 H 6 /Cu h 


410 


240 


0.46 


0.9 


E:Cu ! 


437(1) 


450 






0.17ML Na/Ag 3 


-65* 


-260 






Ar/Ag e 


62(2) 


-1(3) 


0.42(1) 


0.46(4) 


Kr/Ag e 


62(2) 


-08(2) 


0.42(1) 


0.44(4) 


Xe/Ag e 


62(2) 


-52(2) 


0.42(1) 


0.42(6) 


Xe/Ag 1 


67 


-52 


ratio = 


= 1.00(15) 


E:Ag m 


64(1) 


71 







"Ref [108 indicates that the surface state shifts down with Li 
adsorption and disappears (into the bulk) before a half 
monolayer. 6 RefQM] c Ref. H>] d Ref. [UJ e Ref. [TTU] 'Ref. EE]] 
9 Ref. [TT2] ''Ref. [113] 'Ref. EQl^Ref EBl *Ref [HSl ! Ref. UM 
m Ref. Hi 



been recent studies, usi ng t wo-photon photoemission of 
all the alkalis on Cu(llipHand on Ag(lll^SIl they con- 
firm the downward shift but offer little additional quan- 
titative information on the coverage dependence of the 
shifts. To address this shift quantitatively with DFT for 
Na/Cu(lll) wh ile k eeping a manageable cell size, Car- 
avati and TrionP^I U ged a je llium -like model with one- 
dimensional Chulkov potentiate anc j found — e(0) to 
be colinear for coverages = 0, 0.06, and 0.14, of the 
form —0.303 — 2.20, with a slightly smaller negative slope 
when © = 0.25 was also included. For noble-gas adsor- 
bates, the shifts are large while the increase in effective 
mass is small, implying that more than just field effects 
are involved. 

In passing, we discuss a specific complication in deci- 
phering the tabulated numbers: the lateral dipolar in- 
teraction. This effect can influence the field at other 
sites. Most significantly, as alkali atoms get close to 
each other, the direct dipolar repulsion becomes more 
important than the indirect, surface-state- mediated in- 
teraction. This happens when the dipole is large enough 



to produce a significant shift in the surface stateP^ 

We note that the maximum downwards shifts possi- 
ble for metallic surface states are approximatively set 
by the E = value of the minimum surface-state en- 
ergy €f — e E=0 (0), for example, as measured in Ref. [30] 
and reported in Table IV (as the 'clean' entry in the 
'Eisubstrate' rows). The maximum possible upward shift 
in e E ^°(0) is instead approximatively given by the differ- 
ence between value of e E=0 (0) and the energy e bulk (i) = 
90 0meV o f the bulk state at the bottom of the Cu L- 
gap 1 125 * 12 -^ (since the overlap in energies will convert the 
surface state to a surface-state resonance). The maxi- 
mum upward shifts of the surface state for Cu can there- 
fore be estimated by e E=a (0) - e hulk (L) « 460meV. 

When making precise use of the band shifts, one must 
take into account the temperature, as discussed in detail 
for t he ( 111) faces of the three noble metals by Paniago 
et alP3ln particular, they show e F - e E=0 (0) = -(75 ± 
5) meV +(0.17 meV/K) T for Ag(lll). There is a small 
increase also in to, from 0.43 ± 0.04 at 65 K to 0.45 ± 
0.04 at 294 K. For Cu(lll) the increase is comparable, 
with linear coefficient (0.18 ± 0.01) meV/K. 

Overall, it is noteworthy that shifts far larger than 
those reported at the end of Subsection |IVB| are seen. 
The adsorption-induced shifts are indeed larger than the 
shifts of 50 meV which formed the abscissa limits in 
Fig. [5] Hence, the range of field strengths investigated in 
our study are physically sensible. 



C. Decoupling method for adsorbed molecules 

Our study of response to an external field has been 
aided by the decoupling method described and tested in 
Section|TT] This method can also be used to study shift in 
surface-state energy produced by an adsorbed molecule 
or atom. This requires that the two KS states corre- 
sponding to the SL states do not hybridize with other 
bulk or molecular KS states. For systems with inversion 
symmetry in the plane (in terms of the basis vectors), 
like benzene on Cu(lll)J^we only need to decouple the 
states using an 0(2) rotation. 

As a test case, we consider benzene on Cu(lll) in a 
3x3 periodic unit cell and a 6 layer slab. Using DFT with 
the vdW-DF2 functional to capture the non-local corre- 
lation essential to the binding of this system,^ we deter- 
mine a binding separation of 3.5A, and a binding energy 
of 0.48 eV. Except for the use of vdW-DF2, the details of 
this calculation are the same as in our previously reported 
calculation with vdW-DFlP^l After decoupling the two 
KS states, we find a coupling of f2 = 320meV and a de- 
tuning of W = 16meV. Since for the adsorbed-molecule 
case only one of the surfaces is perturbed, the energy 
difference between the decoupled states 2W — 32meV 
equals the shift in surface state energy. 

In this calculation, the decoupling method enabled us 
to extract the shift in surface state energy in a tradi- 
tional slab calculation using far fewer layers than what 
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would be required for the KS to fully decouple because 
of the (weakly) broken symmetry. We have also used the 
method to study the surface shift induced by adsorbed 
anthraquinone chains on Cu(lll) for varying coverages.^ 

VI. CONCLUSIONS 

Using DFT calculations, we have shown that an elec- 
tric field perpendicular to a metal surface with a metallic 
Shockley surface state linearly shifts the bottom of this 
state (relative to the Fermi energy) for physically plau- 
sible field strengths. We have computed the value of the 
linear proportionality constant, as well as that associated 
with the field-induced change in the curvature of the dis- 
persion relation, that is, in the effective mass. 

The decoupling method presented here should be use- 
ful when studying the response in the surface-state dis- 
persion to a perturbation that does not destroy the 
surface-state character. It is, for example, relevant for 
a study of the surface response arising from organic over- 
layers weakly coupled to the surface or dilute overlayers 
of chemisorbed atoms. 

More generally, the MaxLoc analysis could be useful 
for decoupling states which arise at different spatial lo- 



cations and which hybridize under the asumption that an 
infinite time is available to create the coupling. The char- 
acterization of the coupling in such systems can be used 
to calculate tunneling rates and oscillation frequencies. 

It is hard to overemphasize the significance of acquiring 
the ability to control the Fermi wavelength in order to 
manipulate and engineer surface structures determined 
by interactions mediated by surface states. With a strong 
enough field, one could in principle manipulate channels 
in a manner reminiscent of Repp's resonator^^ and so 
dynamically direct an atom flow on surface. 
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